Prettier maps in R and a first step towards interactive maps
Loading some new packages
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✓ ggplot2 3.3.5 ✓ purrr 0.3.4
## ✓ tibble 3.1.3 ✓ dplyr 1.0.7
## ✓ tidyr 1.1.3 ✓ stringr 1.4.0
## ✓ readr 2.0.0 ✓ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(pander)
library(sf)
## Linking to GEOS 3.8.1, GDAL 3.2.1, PROJ 7.2.1
library(terra)
## terra version 1.3.22
##
## Attaching package: 'terra'
## The following object is masked from 'package:pander':
##
## wrap
## The following object is masked from 'package:dplyr':
##
## src
library(units)
## udunits database from /Library/Frameworks/R.framework/Versions/4.1/Resources/library/units/share/udunits/udunits2.xml
library(ggmap)
## Google's Terms of Service: https://cloud.google.com/maps-platform/terms/.
## Please cite ggmap if you use it! See citation("ggmap") for details.
##
## Attaching package: 'ggmap'
## The following object is masked from 'package:terra':
##
## inset
library(cartogram)
##
## Attaching package: 'cartogram'
## The following object is masked from 'package:terra':
##
## cartogram
library(patchwork)
##
## Attaching package: 'patchwork'
## The following object is masked from 'package:terra':
##
## area
library(tmap)
library(viridis)
## Loading required package: viridisLite
library(tigris)
## To enable
## caching of data, set `options(tigris_use_cache = TRUE)` in your R script or .Rprofile.
library(ggspatial)
library(plotly)
##
## Attaching package: 'plotly'
## The following object is masked from 'package:ggmap':
##
## wind
## The following object is masked from 'package:ggplot2':
##
## last_plot
## The following object is masked from 'package:stats':
##
## filter
## The following object is masked from 'package:graphics':
##
## layout
Load the data
#base.dir <- "/Users/mattwilliamson/Google Drive/My Drive/TEACHING/Intro_Spatial_Data_R/Data/session18/"
base.dir <- "/Users/matthewwilliamson/Downloads/session18/"
land.value <- rast(paste0(base.dir,"Regval.tif"))
human.mod <- rast(paste0(base.dir,"hmi.tif"))
mammal.rich <- catalyze(rast(paste0(base.dir,"Mammals_total_richness.tif")))[[2]]
idaho <- states() %>%
filter(., STUSPS == "ID")
##
|
| | 0%
|
|= | 1%
|
|= | 2%
|
|== | 2%
|
|== | 3%
|
|== | 4%
|
|=== | 4%
|
|==== | 6%
|
|===== | 7%
|
|====== | 8%
|
|======= | 10%
|
|======= | 11%
|
|======== | 11%
|
|=========== | 16%
|
|============ | 16%
|
|================ | 23%
|
|====================== | 31%
|
|====================== | 32%
|
|========================= | 36%
|
|========================== | 36%
|
|========================== | 37%
|
|========================== | 38%
|
|=========================== | 38%
|
|=========================== | 39%
|
|=============================== | 45%
|
|================================ | 45%
|
|================================ | 46%
|
|====================================== | 54%
|
|====================================== | 55%
|
|======================================== | 57%
|
|======================================== | 58%
|
|========================================= | 58%
|
|========================================= | 59%
|
|========================================== | 59%
|
|========================================== | 60%
|
|============================================== | 65%
|
|============================================== | 66%
|
|=============================================== | 67%
|
|=============================================== | 68%
|
|================================================ | 68%
|
|================================================ | 69%
|
|================================================= | 69%
|
|================================================= | 70%
|
|================================================== | 72%
|
|=================================================== | 72%
|
|=================================================== | 73%
|
|==================================================== | 74%
|
|==================================================== | 75%
|
|===================================================== | 75%
|
|====================================================== | 77%
|
|====================================================== | 78%
|
|======================================================= | 78%
|
|======================================================= | 79%
|
|======================================================== | 80%
|
|=========================================================== | 84%
|
|=========================================================== | 85%
|
|============================================================ | 85%
|
|================================================================ | 91%
|
|================================================================ | 92%
|
|===================================================================== | 98%
|
|===================================================================== | 99%
|
|======================================================================| 99%
|
|======================================================================| 100%
land.value.crop <- terra::crop(land.value, terra::project(vect(idaho), crs(land.value)))
human.mod.crop <- terra::crop(human.mod, terra::project(vect(idaho), crs(human.mod)))
mammal.rich.crop <- terra::crop(mammal.rich, terra::project(vect(idaho), crs(mammal.rich)))
land.value.project <- terra::resample(
terra::project(land.value.crop, crs(human.mod.crop)),
human.mod.crop)
mammal.rich.project <- terra::resample(
terra::project(mammal.rich.crop, crs(human.mod.crop)),
human.mod.crop)
rast.stack <- c(human.mod.crop, land.value.project, mammal.rich.project)
idaho <- idaho %>%
st_transform(., crs = st_crs(rast.stack))
idaho.counties <- counties(state = "ID") %>%
st_transform(., crs = st_crs(rast.stack))
##
|
| | 0%
|
|= | 1%
|
|= | 2%
|
|== | 3%
|
|== | 4%
|
|=== | 4%
|
|=== | 5%
|
|==== | 5%
|
|==== | 6%
|
|===== | 6%
|
|===== | 7%
|
|===== | 8%
|
|====== | 8%
|
|====== | 9%
|
|======= | 9%
|
|======= | 10%
|
|======= | 11%
|
|======== | 11%
|
|======== | 12%
|
|========= | 12%
|
|========= | 13%
|
|========= | 14%
|
|========== | 14%
|
|========== | 15%
|
|=========== | 15%
|
|=========== | 16%
|
|============ | 17%
|
|============= | 18%
|
|============= | 19%
|
|============== | 19%
|
|============== | 20%
|
|============== | 21%
|
|=============== | 21%
|
|=============== | 22%
|
|================ | 22%
|
|================ | 23%
|
|================ | 24%
|
|================= | 24%
|
|================= | 25%
|
|================== | 25%
|
|================== | 26%
|
|=================== | 26%
|
|=================== | 27%
|
|=================== | 28%
|
|==================== | 28%
|
|==================== | 29%
|
|===================== | 30%
|
|===================== | 31%
|
|====================== | 31%
|
|====================== | 32%
|
|======================= | 32%
|
|======================= | 33%
|
|======================== | 34%
|
|======================== | 35%
|
|========================= | 35%
|
|========================= | 36%
|
|========================== | 36%
|
|========================== | 37%
|
|========================== | 38%
|
|=========================== | 38%
|
|=========================== | 39%
|
|============================ | 39%
|
|============================ | 40%
|
|============================ | 41%
|
|============================= | 41%
|
|============================= | 42%
|
|============================== | 42%
|
|============================== | 43%
|
|=============================== | 44%
|
|================================ | 45%
|
|================================ | 46%
|
|================================= | 47%
|
|================================= | 48%
|
|================================== | 48%
|
|================================== | 49%
|
|=================================== | 50%
|
|=================================== | 51%
|
|==================================== | 51%
|
|==================================== | 52%
|
|===================================== | 52%
|
|===================================== | 53%
|
|====================================== | 54%
|
|======================================= | 55%
|
|======================================= | 56%
|
|======================================== | 57%
|
|========================================= | 58%
|
|========================================= | 59%
|
|========================================== | 59%
|
|========================================== | 60%
|
|=========================================== | 61%
|
|=========================================== | 62%
|
|============================================ | 62%
|
|============================================ | 63%
|
|============================================= | 64%
|
|============================================= | 65%
|
|============================================== | 65%
|
|============================================== | 66%
|
|=============================================== | 67%
|
|=============================================== | 68%
|
|================================================ | 68%
|
|================================================ | 69%
|
|================================================= | 70%
|
|================================================= | 71%
|
|================================================== | 71%
|
|================================================== | 72%
|
|=================================================== | 72%
|
|=================================================== | 73%
|
|==================================================== | 74%
|
|===================================================== | 75%
|
|===================================================== | 76%
|
|====================================================== | 77%
|
|======================================================= | 78%
|
|======================================================= | 79%
|
|======================================================== | 79%
|
|======================================================== | 80%
|
|======================================================== | 81%
|
|========================================================= | 81%
|
|========================================================= | 82%
|
|========================================================== | 82%
|
|========================================================== | 83%
|
|=========================================================== | 84%
|
|=========================================================== | 85%
|
|============================================================ | 85%
|
|============================================================ | 86%
|
|============================================================= | 86%
|
|============================================================= | 87%
|
|============================================================== | 88%
|
|============================================================== | 89%
|
|=============================================================== | 89%
|
|=============================================================== | 90%
|
|=============================================================== | 91%
|
|================================================================ | 91%
|
|================================================================ | 92%
|
|================================================================= | 92%
|
|================================================================= | 93%
|
|================================================================= | 94%
|
|================================================================== | 94%
|
|================================================================== | 95%
|
|=================================================================== | 95%
|
|=================================================================== | 96%
|
|==================================================================== | 97%
|
|===================================================================== | 98%
|
|===================================================================== | 99%
|
|======================================================================| 99%
|
|======================================================================| 100%
idaho.pas <- st_read(paste0(base.dir,"westernPAs.shp")) %>%
filter(., Stat_Nm == "ID" ) %>%
st_transform(., crs = st_crs(rast.stack))
## Reading layer `westernPAs' from data source
## `/Users/matthewwilliamson/Downloads/session18/westernPAs.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 8341 features and 34 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -2395532 ymin: 996849.8 xmax: -512593.7 ymax: 3172010
## Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic_USGS_version
id.census <- tidycensus:: get_acs(geography = "county",
variables = c(medianincome = "B19013_001",
pop = "B01003_001"),
state = c("ID"),
year = 2018,
key = key,
geometry = TRUE) %>%
st_transform(., crs = st_crs(rast.stack)) %>%
dplyr::select(-moe) %>%
spread(variable, estimate)
## Getting data from the 2014-2018 5-year ACS
## Downloading feature geometry from the Census website. To cache shapefiles for use in future sessions, set `options(tigris_use_cache = TRUE)`.
##
|
| | 0%
|
|=== | 5%
|
|==== | 5%
|
|==== | 6%
|
|======== | 12%
|
|========= | 12%
|
|========= | 13%
|
|============= | 19%
|
|============== | 19%
|
|============== | 20%
|
|=============== | 21%
|
|================ | 23%
|
|================ | 24%
|
|================= | 24%
|
|================= | 25%
|
|===================== | 29%
|
|===================== | 30%
|
|===================== | 31%
|
|====================== | 31%
|
|========================== | 37%
|
|========================== | 38%
|
|=========================== | 38%
|
|=============================== | 44%
|
|=============================== | 45%
|
|================================ | 45%
|
|=================================== | 51%
|
|==================================== | 51%
|
|==================================== | 52%
|
|===================================== | 52%
|
|======================================== | 58%
|
|========================================= | 58%
|
|=========================================== | 62%
|
|============================================ | 62%
|
|============================================ | 63%
|
|============================================== | 65%
|
|================================================ | 69%
|
|================================================== | 72%
|
|=================================================== | 72%
|
|=================================================== | 73%
|
|======================================================= | 79%
|
|======================================================== | 79%
|
|======================================================== | 80%
|
|============================================================ | 85%
|
|============================================================ | 86%
|
|================================================================ | 91%
|
|================================================================ | 92%
|
|================================================================= | 93%
|
|===================================================================== | 99%
|
|======================================================================| 100%
#idaho.elevation <- get_elev_raster(locations = idaho, z = 8, expand = 10000)
Building a map with a few layers
bg <- ggmap::get_map(as.vector(st_bbox(st_transform(idaho, 4326))), zoom = 7)
## Source : http://tile.stamen.com/terrain/7/22/43.png
## Source : http://tile.stamen.com/terrain/7/23/43.png
## Source : http://tile.stamen.com/terrain/7/24/43.png
## Source : http://tile.stamen.com/terrain/7/22/44.png
## Source : http://tile.stamen.com/terrain/7/23/44.png
## Source : http://tile.stamen.com/terrain/7/24/44.png
## Source : http://tile.stamen.com/terrain/7/22/45.png
## Source : http://tile.stamen.com/terrain/7/23/45.png
## Source : http://tile.stamen.com/terrain/7/24/45.png
## Source : http://tile.stamen.com/terrain/7/22/46.png
## Source : http://tile.stamen.com/terrain/7/23/46.png
## Source : http://tile.stamen.com/terrain/7/24/46.png
## Source : http://tile.stamen.com/terrain/7/22/47.png
## Source : http://tile.stamen.com/terrain/7/23/47.png
## Source : http://tile.stamen.com/terrain/7/24/47.png
ggmap(bg) +
geom_sf(data = st_transform(id.census, 4326), mapping = aes(fill = medianincome), inherit.aes = FALSE) +
geom_sf(data=st_transform(idaho.pas, 4326), color="yellow", fill=NA, inherit.aes = FALSE) +
scale_fill_continuous() +
coord_sf(crs = st_crs(4326))
## Coordinate system already present. Adding new coordinate system, which will replace the existing one.

ggmap(bg) +
geom_sf(data = st_transform(id.census, 4326), mapping = aes(fill = medianincome), inherit.aes = FALSE) +
geom_sf(data=st_transform(idaho.pas, 4326), color="yellow", fill=NA, inherit.aes = FALSE) +
scale_fill_continuous() +
coord_sf(crs = st_crs(4326)) +
theme(legend.direction = "horizontal", legend.position = "bottom", legend.justification = "center") +
theme_bw()
## Coordinate system already present. Adding new coordinate system, which will replace the existing one.

ggmap(bg) +
geom_sf(data = st_transform(id.census, 4326), mapping = aes(fill = medianincome), inherit.aes = FALSE) +
geom_sf(data=st_transform(idaho.pas, 4326), color="yellow", fill=NA, inherit.aes = FALSE) +
scale_fill_continuous() +
coord_sf(crs = st_crs(4326)) +
annotation_scale(location = "tl") +
annotation_north_arrow(location = "br", which_north = "true") +
ggtitle("PAs in Idaho") +
theme(legend.direction = "horizontal", legend.position = "bottom", legend.justification = "center") +
theme_bw()
## Coordinate system already present. Adding new coordinate system, which will replace the existing one.
## Scale on map varies by more than 10%, scale bar may be inaccurate

main.map <- ggmap(bg) +
geom_sf(data = st_transform(id.census, 4326), mapping = aes(fill = medianincome), inherit.aes = FALSE) +
geom_sf(data=st_transform(idaho.pas, 4326), color="yellow", fill=NA, inherit.aes = FALSE) +
scale_fill_continuous() +
coord_sf(crs = st_crs(4326)) +
annotation_scale(location = "tl") +
annotation_north_arrow(location = "br", which_north = "true") +
ggtitle("PAs in Idaho") +
theme(legend.direction = "horizontal", legend.position = "bottom", legend.justification = "center") +
theme_bw()
not.conus <- c("AK", "HI", "DC", "MP", "GU", "VI", "AS", "PR")
conus <- states() %>%
filter(., !(STUSPS %in% not.conus)) %>%
st_transform(., 4326)
bbox <- st_as_sfc(st_bbox(st_transform(idaho, 4326)))
inset.map <- ggplot(conus)+
geom_sf(fill="lightgray", color="black") +
geom_sf(data = st_as_sfc(st_bbox(st_transform(idaho, 4326))),fill=NA, color = "red") +
theme_blank()
complete <- main.map + inset_element(inset.map, left = 0.6, bottom = 0.6, right = 1, top = 1, align_to = "full")
ggsave("insetmap.png", plot=complete)

p <- ggplot() +
geom_sf(data = st_transform(id.census, 4326), mapping = aes(fill = medianincome), inherit.aes = FALSE) +
geom_sf(data=st_transform(idaho.pas, 4326), color="yellow", fill=NA, inherit.aes = FALSE) +
scale_fill_continuous()
ggplotly(p)